1 Data and package

library(scran)
library(SingleCellExperiment)
library(scater)
library(scattermore)
library(moon)
library(ggplot2)
library(ggthemes)
library(ggpubr)
library(reshape2)
library(dplyr)
library(stringr)
library(pheatmap)
library(CellChat)
library(RColorBrewer)
library(gridExtra)
coldata_wilks <- readRDS("results/wilk_results/coldata_wilks.rds")
meta_wilks <- readRDS("results/wilk_results/meta_wilks.rds")
# CCI results (Cellchat)
cellchat_res_list <- readRDS("results/wilk_results/wilk_cellchat_res_list.rds")
rankNet_byCellType <- function(object, slot.name = "netP", 
                               x.rotation = 90, title = NULL, color.use = NULL, 
                               bar.w = 0.75, font.size = 8) 
{
  object1 <- methods::slot(object, slot.name)
  prob1 = object1$prob
  df <- melt(apply(prob1, 3, function(x) {
    df <- melt(x)
    colnames(df) <- c("Ligand", "Receptor", "value")
    df
  }))
  df <- df[, c("Ligand", "Receptor", "L1", "value")]
  colnames(df)[3] <- "Pathway"
  return(df)
  
  
}

2 Overall pattern analysis

rankNet_byCellType_list <- lapply(cellchat_res_list, rankNet_byCellType)

rankNet_byCellType_list <- melt(rankNet_byCellType_list)
rankNet_byCellType_list$Ligand_group <- unlist(lapply(strsplit(as.character(rankNet_byCellType_list$Ligand), 
                                                               "_"), "[[", 1))
rankNet_byCellType_list$Receptor_group <- unlist(lapply(strsplit(as.character(rankNet_byCellType_list$Receptor), 
                                                                 "_"), "[[", 1))


saveRDS(rankNet_byCellType_list, file = "results/wilk_results/rankNet_byCellType_list_wilks.rds")
severity_color <- c("#2ca02c", "#FFD92F", "#7570B3")
names(severity_color) <- c("Healthy", "Moderate", "Severe")
rankNet_byGroup_agg <- aggregate(rankNet_byCellType_list$value, 
                                 list(rankNet_byCellType_list$Ligand_group,
                                      rankNet_byCellType_list$Receptor_group,
                                      rankNet_byCellType_list$L1,
                                      rankNet_byCellType_list$Pathway),
                                 sum)


colnames(rankNet_byGroup_agg) <- c("Ligand_group", 
                                   "Receptor_group",
                                   "sample",
                                   "Pathway",
                                   "value")
features <- paste(rankNet_byGroup_agg$Ligand_group,
                  rankNet_byGroup_agg$Receptor_group,
                  rankNet_byGroup_agg$Pathway, sep = "_")

rankNet_byGroup_agg$features <- features
rankNet_byGroup_agg_all <- dcast2(rankNet_byGroup_agg, 
                                  features ~ sample, 
                                  fun.aggregate = sum, value.var = "value")
rankNet_byGroup_agg_all <- rankNet_byGroup_agg_all[rowSums(rankNet_byGroup_agg_all) > 0, ]
rankNet_byGroup_agg_all <- rankNet_byGroup_agg_all[rowSums(rankNet_byGroup_agg_all!=0) > 2, ]

2.1 Feature selction: kruskal test

kruskal_pvalue <- list()
for (i in 1:nrow(rankNet_byGroup_agg_all)) {
   # if (i %% 100 == 0) print(i)
    
    kruskal_res <- try(kruskal.test(unlist(rankNet_byGroup_agg_all[i,]) ~ meta_wilks[colnames(rankNet_byGroup_agg_all), ]$Condition), silent = TRUE)
    kruskal_pvalue[[i]] <- try(kruskal_res$p.value, silent = TRUE)
    
}

kruskal_pvalue <- lapply(kruskal_pvalue, function(x) {
    if (class(x) == "try-error") {
        x <- NULL
    }
    x
})
names(kruskal_pvalue) <- rownames(rankNet_byGroup_agg_all)
kruskal_pvalue <- unlist(kruskal_pvalue)

kruskal_pvalue <- p.adjust(kruskal_pvalue, method = "BH")
saveRDS(kruskal_pvalue, file = "results/wilk_results/CCI_kruskal_pvalue_condition_wilk.rds")

2.2 PCA

pca_patient <- prcomp(t(-1/log(rankNet_byGroup_agg_all[names(kruskal_pvalue[kruskal_pvalue < 0.2]),])), 
                      scale. = TRUE, center = TRUE)
library(ggrepel)
pca1 <- ggplot(data.frame(pca_patient$x), aes(x = pca_patient$x[, 1],
                                              y = pca_patient$x[, 2],
                                              color = meta_wilks[rownames(pca_patient$x),]$Condition)) +
    geom_point(size = 4, alpha = 0.8) +
    # geom_text_repel(aes(label = rownames(pca_patient$x))) +
    theme_yx() +
    theme(aspect.ratio = 1) +
    scale_color_manual(values = severity_color) +
    xlab("PCA1") +
    ylab("PCA2") +
    labs(color = "")

pca2 <- ggplot(data.frame(pca_patient$x), aes(x = pca_patient$x[, 1],
                                              y = pca_patient$x[, 3],
                                              color = meta_wilks[rownames(pca_patient$x),]$Condition)) +
    geom_point(size = 3, alpha = 0.8) +
    # geom_text_repel(aes(label = rownames(pca_patient$x))) +
    theme_yx() +
    theme(aspect.ratio = 1) +
    scale_color_manual(values = severity_color) +
    xlab("PCA1") +
    ylab("PCA3") +
    labs(color = "")

pca3 <- ggplot(data.frame(pca_patient$x), aes(x = pca_patient$x[, 2],
                                              y = pca_patient$x[, 3],
                                              color = meta_wilks[rownames(pca_patient$x),]$Condition)) +
    geom_point(size = 3, alpha = 0.8) +
    # geom_text_repel(aes(label = rownames(pca_patient$x))) +
    theme_yx() +
    theme(aspect.ratio = 1) +
    scale_color_manual(values = severity_color) +
    xlab("PCA2") +
    ylab("PCA3") +
    labs(color = "")

ggarrange(pca1, pca2, pca3, align = "hv", 
          common.legend = TRUE, ncol = 2, nrow = 2)

2.3 Aggregation by samples

aff_mat_bySample <- lapply(split(rankNet_byGroup_agg, rankNet_byGroup_agg$sample),
                           function(x) dcast2(x, Ligand_group~Receptor_group,
                                              fun.aggregate = mean, value.var = "value"))
all_cellTypes <- names(table(coldata_wilks$cell.type.fine))

aff_mat_bySample <- lapply(aff_mat_bySample, function(x) {
    mat <- matrix(0, ncol = length(all_cellTypes), nrow = length(all_cellTypes))
    colnames(mat) <- rownames(mat) <- all_cellTypes
    mat[rownames(x), colnames(x)] <- as.matrix(x)
    mat
})


selected_cellTypes <- all_cellTypes[!all_cellTypes %in% c("intermediate", "unassigned",
                                                         "RBC", "Platelet", "SC & Eosinophil")]
aff_mat_bySample <- lapply(aff_mat_bySample, function(x) {
  x <- x[selected_cellTypes, selected_cellTypes]
  x
})


aff_mat_bySample <- lapply(aff_mat_bySample, function(x) {
    (x - min(x))/(max(x) - min(x))
})


p <- lapply(1:length(aff_mat_bySample), function(i) {
    pheatmap(aff_mat_bySample[[i]],
             cluster_cols = FALSE, 
             cluster_rows = FALSE,
             main = names(aff_mat_bySample)[i],
             color =  colorRampPalette(c("white", 
                                         brewer.pal(n = 7, 
                                                    name = "Reds")))(100))
})

pdf("figures/WilkEtAl/cellchat_CCI_network_sample_byCellType.pdf", 
    width = 20, height = 20)
do.call(grid.arrange, list(grobs = lapply(p, function(x) x$gtable), ncol = 3))
dev.off()
## quartz_off_screen 
##                 2
# 


severe_patients <- rownames(meta_wilks)[meta_wilks$Condition == "Severe"]
aff_mat_severe <- Reduce("+", aff_mat_bySample[names(aff_mat_bySample) %in% severe_patients])/length(severe_patients)

moderate_patients <- rownames(meta_wilks)[meta_wilks$Condition == "Moderate"]
aff_mat_moderate <- Reduce("+", aff_mat_bySample[names(aff_mat_bySample) %in% moderate_patients])/length(moderate_patients)

control_patients <- rownames(meta_wilks)[meta_wilks$Condition == "Healthy"]
aff_mat_control <- Reduce("+", aff_mat_bySample[names(aff_mat_bySample) %in% control_patients])/length(control_patients)



p_severe <- pheatmap(aff_mat_severe, cluster_cols = FALSE, 
                     cluster_rows = FALSE,
                     main = "severe (average across samples)",
                     color =  colorRampPalette(c("white", 
                                                 brewer.pal(n = 7, 
                                                            name = "Reds")))(100),
                     breaks = seq(0, max(aff_mat_severe), max(aff_mat_severe)/100))

library(RColorBrewer)
p_moderate <- pheatmap(aff_mat_moderate, 
                       cluster_cols = FALSE, 
                       cluster_rows = FALSE,
                       main = "moderate (average across samples)",
                       color =  colorRampPalette(c("white", 
                                                   brewer.pal(n = 7, 
                                                              name = "Reds")))(100),
                       breaks = seq(0, max(aff_mat_moderate), max(aff_mat_moderate)/100))

p_control <- pheatmap(aff_mat_control, 
                      cluster_cols = FALSE, 
                      cluster_rows = FALSE,
                      main = "control (average across samples)",
                      color =  colorRampPalette(c("white", 
                                                  brewer.pal(n = 7, 
                                                             name = "Reds")))(100),
                      breaks = seq(0, max(aff_mat_control), max(aff_mat_control)/100))

pdf("figures/WilkEtAl/cellchat_CCI_network_byCondition_noScale.pdf", 
    width = 18, height = 6)
do.call(grid.arrange, list(grobs = list(p_control$gtable,
                                        p_moderate$gtable,
                                        p_severe$gtable), ncol = 3))
dev.off()
## quartz_off_screen 
##                 2
aff_mat_diff <- aff_mat_severe - aff_mat_moderate

pheatmap(aff_mat_diff,
         cluster_cols = FALSE, 
         cluster_rows = FALSE,
         color =  colorRampPalette(c("blue", "white", "red"))(100),
         breaks = c(seq(min(aff_mat_diff), 0, (0 - min(aff_mat_diff))/50),
                    seq(0.01, max(aff_mat_diff), (max(aff_mat_diff))/50)),
         main = "server - moderate (wilks et al.)",
         #file = "figures/WilkEtAl/cellchat_CCI_network_byCondition_diff_severe.pdf",
         width = 8,
         height = 7)

aff_mat_diff <- aff_mat_moderate - aff_mat_control


pheatmap(aff_mat_diff,
         cluster_cols = FALSE, 
         cluster_rows = FALSE,
         color =  colorRampPalette(c("blue", "white", "red"))(100),
         breaks = c(seq(min(aff_mat_diff), 0, (0 - min(aff_mat_diff))/50),
                    seq(0.01, max(aff_mat_diff), (max(aff_mat_diff))/50)),
         main = "moderate - control (wilks et al.)",
         #file = "figures/WilkEtAl/cellchat_CCI_network_byCondition_diff_moderate.pdf",
         width = 8,
         height = 7)

aff_mat_diff <- aff_mat_severe - aff_mat_control


pheatmap(aff_mat_diff,
         cluster_cols = FALSE, 
         cluster_rows = FALSE,
         color =  colorRampPalette(c("blue", "white", "red"))(100),
         breaks = c(seq(min(aff_mat_diff), 0, (0 - min(aff_mat_diff))/50),
                    seq(0.01, max(aff_mat_diff), (max(aff_mat_diff))/50)),
         main = "severe - control (wilks et al.)",
         #file = "figures/WilkEtAl/cellchat_CCI_network_byCondition_diff_severe_control.pdf",
         width = 8,
         height = 7)

cellType_col <- moon_pal(22)[-c(1:2)][c(20, 1:19)]
cellType_col <- append(cellType_col, c("grey80", "grey30"))
names(cellType_col) <- c(names(table(coldata_wilks$cell.type)), c("intermediate", "unassigned"))
saveRDS(cellType_col, file = "results/cellType_col")
mat <- aff_mat_control


cci_control <- melt(as.matrix(mat))
colnames(cci_control) <- c("Ligand", "Receptor", "n")
library(igraph)
library(ggraph)
g <- graph_from_data_frame(data.frame(cci_control))
E(g)$weights <- ifelse(cci_control$n == 0,
                       1e-10, abs(cci_control$n))

V(g)$color <- cellType_col[V(g)$name]
# pdf("figures/WilkEtAl/cellchat_CCI_network_byCondition_Control_network.pdf",
#     width = 8,
#     height = 6)
plot(g, 
     vertex.size = 20,
     vertex.color = V(g)$color,
     vertex.label.color = "black",
     vertex.label.cex = 1,
     edge.width = E(g)$weights * 20,
     edge.arrow.size = log(1/E(g)$weights)/80,
     #edge.arrow.size = 0.01,
     #edge.color = E(g)$sign,
     edge.curved = 0.3,
     layout = layout_in_circle,
     main = "Control")

##   [1] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [19] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [37] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [55] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [73] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [91] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [109] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [127] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [145] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [163] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [181] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [199] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [217] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [235] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [253] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [271] 0.3 0.3
# dev.off()
mat <- aff_mat_moderate


cci_moderate <- melt(as.matrix(mat))
colnames(cci_moderate) <- c("Ligand", "Receptor", "n")
library(igraph)
library(ggraph)
g <- graph_from_data_frame(data.frame(cci_moderate))
E(g)$weights <- ifelse(cci_moderate$n == 0,
                       1e-10, abs(cci_moderate$n))
#E(g)$sign <- ifelse(sign(cci_moderate$n) == 1, "#d62728", "#1f77b4")
g <- delete_edges(g, E(g)[cci_moderate$n == 0])

V(g)$color <- cellType_col[V(g)$name]
# pdf("figures/WilkEtAl/cellchat_CCI_network_byCondition_moderate_network.pdf",
#     width = 8,
#     height = 6)
plot(g, 
     vertex.size = 20,
     vertex.color = V(g)$color,
     vertex.label.color = "black",
     vertex.label.cex = 1,
     edge.width = E(g)$weights * 20,
     edge.arrow.size = log(1/E(g)$weights)/20,
     #edge.arrow.size = 0.01,
     #edge.color = E(g)$sign,
     edge.curved = 0.3,
     layout = layout_in_circle,
     main = "moderate")

##   [1] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [19] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [37] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [55] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [73] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [91] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [109] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [127] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [145] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
# dev.off()
mat <- aff_mat_severe

cci_severe <- melt(as.matrix(mat))
colnames(cci_severe) <- c("Ligand", "Receptor", "n")
library(igraph)
library(ggraph)
g <- graph_from_data_frame(data.frame(cci_severe))

E(g)$weights <- ifelse(cci_severe$n == 0,
                       1e-10, abs(cci_severe$n))

g <- delete_edges(g, E(g)[cci_severe$n == 0])
V(g)$color <- cellType_col[V(g)$name]
# pdf("figures/WilkEtAl/cellchat_CCI_network_byCondition_severe_network.pdf",
#     width = 8,
#     height = 6)
plot(g, 
     vertex.size = 20,
     vertex.color = V(g)$color,
     vertex.label.color = "black",
     vertex.label.cex = 1,
     edge.width = E(g)$weights * 20,
     edge.arrow.size = log(1/E(g)$weights)/20,
     #edge.arrow.size = 0.01,
     #edge.color = E(g)$sign,
     edge.curved = 0.3,
     layout = layout_in_circle,
     main = "severe")

##   [1] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [19] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [37] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [55] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [73] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
##  [91] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [109] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [127] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [145] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [163] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [181] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [199] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [217] 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3
## [235] 0.3 0.3 0.3 0.3
# dev.off()

3 Session Info

sessionInfo()
## R version 4.0.2 RC (2020-06-20 r78727)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Catalina 10.15.7
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRblas.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib
## 
## locale:
## [1] en_AU.UTF-8/en_AU.UTF-8/en_AU.UTF-8/C/en_AU.UTF-8/en_AU.UTF-8
## 
## attached base packages:
## [1] parallel  stats4    stats     graphics  grDevices utils     datasets 
## [8] methods   base     
## 
## other attached packages:
##  [1] ggraph_2.0.3                igraph_1.1.0               
##  [3] ggrepel_0.8.2               gridExtra_2.3              
##  [5] RColorBrewer_1.1-2          CellChat_0.0.1             
##  [7] bigmemory_4.5.36            pheatmap_1.0.12            
##  [9] stringr_1.4.0               dplyr_1.0.2                
## [11] reshape2_1.4.4              ggpubr_0.3.0               
## [13] ggthemes_4.2.0              moon_0.1.0                 
## [15] scattermore_0.6             scater_1.16.1              
## [17] ggplot2_3.3.2               scran_1.16.0               
## [19] SingleCellExperiment_1.10.1 SummarizedExperiment_1.18.1
## [21] DelayedArray_0.14.0         matrixStats_0.56.0         
## [23] Biobase_2.48.0              GenomicRanges_1.40.0       
## [25] GenomeInfoDb_1.24.2         IRanges_2.22.2             
## [27] S4Vectors_0.26.1            BiocGenerics_0.34.0        
## 
## loaded via a namespace (and not attached):
##   [1] readxl_1.3.1              backports_1.1.8          
##   [3] circlize_0.4.10           systemfonts_0.2.3        
##   [5] NMF_0.30.1                plyr_1.8.6               
##   [7] BiocParallel_1.22.0       listenv_0.8.0            
##   [9] gridBase_0.4-7            digest_0.6.25            
##  [11] foreach_1.5.0             htmltools_0.5.0          
##  [13] viridis_0.5.1             ggalluvial_0.12.0        
##  [15] magrittr_1.5              cluster_2.1.0            
##  [17] doParallel_1.0.15         openxlsx_4.1.5           
##  [19] limma_3.44.3              graphlayouts_0.7.0       
##  [21] sna_2.5                   ComplexHeatmap_2.4.2     
##  [23] globals_0.12.5            svglite_1.2.3.2          
##  [25] colorspace_1.4-1          haven_2.3.1              
##  [27] xfun_0.18                 crayon_1.3.4             
##  [29] RCurl_1.98-1.2            jsonlite_1.6.1           
##  [31] bigmemory.sri_0.1.3       iterators_1.0.12         
##  [33] glue_1.4.1                polyclip_1.10-0          
##  [35] pals_1.6                  registry_0.5-1           
##  [37] gtable_0.3.0              zlibbioc_1.34.0          
##  [39] XVector_0.28.0            GetoptLong_1.0.0         
##  [41] car_3.0-8                 BiocSingular_1.4.0       
##  [43] future.apply_1.5.0        shape_1.4.4              
##  [45] maps_3.3.0                abind_1.4-5              
##  [47] scales_1.1.1              edgeR_3.30.3             
##  [49] rngtools_1.5              bibtex_0.4.2.2           
##  [51] rstatix_0.6.0             Rcpp_1.0.4.6             
##  [53] viridisLite_0.3.0         xtable_1.8-4             
##  [55] clue_0.3-57               reticulate_1.16          
##  [57] dqrng_0.2.1               mapproj_1.2.7            
##  [59] foreign_0.8-80            rsvd_1.0.3               
##  [61] FNN_1.1.3                 ellipsis_0.3.1           
##  [63] farver_2.0.3              pkgconfig_2.0.3          
##  [65] locfit_1.5-9.4            labeling_0.3             
##  [67] tidyselect_1.1.0          rlang_0.4.9              
##  [69] munsell_0.5.0             cellranger_1.1.0         
##  [71] tools_4.0.2               generics_0.0.2           
##  [73] statnet.common_4.3.0      broom_0.7.2              
##  [75] evaluate_0.14             yaml_2.2.1               
##  [77] knitr_1.30                tidygraph_1.2.0          
##  [79] zip_2.0.4                 purrr_0.3.4              
##  [81] dendextend_1.13.4         pbapply_1.4-2            
##  [83] future_1.17.0             compiler_4.0.2           
##  [85] beeswarm_0.2.3            curl_4.3                 
##  [87] png_0.1-7                 ggsignif_0.6.0           
##  [89] tweenr_1.0.1              tibble_3.0.4             
##  [91] statmod_1.4.34            stringi_1.4.6            
##  [93] RSpectra_0.16-0           gdtools_0.2.2            
##  [95] forcats_0.5.0             lattice_0.20-41          
##  [97] Matrix_1.2-18             vctrs_0.3.5              
##  [99] pillar_1.4.4              lifecycle_0.2.0          
## [101] GlobalOptions_0.1.2       BiocNeighbors_1.6.0      
## [103] data.table_1.12.8         cowplot_1.0.0            
## [105] bitops_1.0-6              irlba_2.3.3              
## [107] R6_2.4.1                  network_1.16.0           
## [109] rio_0.5.16                vipor_0.4.5              
## [111] codetools_0.2-16          dichromat_2.0-0          
## [113] MASS_7.3-51.6             assertthat_0.2.1         
## [115] pkgmaker_0.31.1           rjson_0.2.20             
## [117] withr_2.2.0               GenomeInfoDbData_1.2.3   
## [119] hms_0.5.3                 grid_4.0.2               
## [121] coda_0.19-3               tidyr_1.1.2              
## [123] rmarkdown_2.4             DelayedMatrixStats_1.10.0
## [125] carData_3.0-4             ggforce_0.3.1            
## [127] ggbeeswarm_0.6.0